clear all
clf
close all
% =====================================================================
% ====  "Risk-Sharing in Village Economies Revisited", by Bold and Broer

% ====   This programme draws the graph of asymmetry moments as function of
% village size

% ====    Written by Tessa Bold and Tobias Broer

% ====    This version Sep 2021
% =====================================================================
write=0; % set to 1 to write figures to file

currentfolder='C:\Users\tbroer\Dropbox\Research\Current_Projects\tessa\JEEA Replication\Programmes\Main Results and Robustness';
currentfolder1='C:\Users\tbroer\Dropbox\Research\Current_Projects\tessa\JEEA Replication\Programmes\Main Results and Robustness';
%  ================================
%  Figure 4
%  ================================
load([currentfolder '\BuildupMomentsest101010'])

plotasym=1; % set to 1 if you want the asymmetry
date='5_June_21';
delta_calibvec=[0.86,0.91];
rho_calib=1;

count=0;
beta_dcdy_data=0.21 % regression coefficient of dc on dy for Kanzara
beta_relvardcdy_data=0.30 % regression coefficient of dc on dy for Kanzara

figure2=figure
set(figure2,'defaultAxesColorOrder',[[0 0 0];[0 0 0]])
for delta_calib=delta_calibvec
    count=count+1;
[index1 index2] = find(abs(squeeze(Momentsest(2,1,:,:,end))-delta_calib )<0.000000001 & squeeze(Momentsest(3,1,:,:,end))==rho_calib );
vss_set=find(~isnan([squeeze(Momentsest(1,1,index1,index2,1:end))]).*squeeze(Momentsest(1,1,index1,index2,1:end))>0);
vss_max=vss_set(end);
xx=vss_set;
buildupasym=nan(3,length(vss_set),3);

%unconditional
buildupasym(1,:,1)=[squeeze(Momentsest(7,1,index1,index2,[vss_set]))];%./squeeze(Momentsest(21,1,index1,index2,[vss_set]))
buildupasym(2,:,1)=[squeeze(Momentsest(13,1,index1,index2,[vss_set]))-squeeze(Momentsest(14,1,index1,index2,[vss_set]))];
buildupasym(3,:,1)=[(squeeze(Momentsest(13,1,index1,index2,[vss_set]))-squeeze(Momentsest(14,1,index1,index2,[vss_set])))./squeeze(Momentsest(12,1,index1,index2,[vss_set]))];
buildupasym(4,:,1)=[squeeze(Momentsest(6,1,index1,index2,[vss_set]))];
%conditional on village income
buildupasym(1,:,2)=[squeeze(Momentsest(9,1,index1,index2,[vss_set]))];
buildupasym(2,:,2)=[squeeze(Momentsest(16,1,index1,index2,[vss_set]))-squeeze(Momentsest(17,1,index1,index2,[vss_set]))];
buildupasym(3,:,2)=[(squeeze(Momentsest(16,1,index1,index2,[vss_set]))-squeeze(Momentsest(17,1,index1,index2,[vss_set])))];
buildupasym(4,:,2)=[squeeze(Momentsest(15,1,index1,index2,[vss_set]))];
%conditional on group income
buildupasym(1,:,3)=[squeeze(Momentsest(11,1,index1,index2,[vss_set]))];
buildupasym(2,:,3)=[squeeze(Momentsest(19,1,index1,index2,[vss_set]))-squeeze(Momentsest(20,1,index1,index2,[vss_set]))];
buildupasym(3,:,3)=[(squeeze(Momentsest(19,1,index1,index2,[vss_set]))-squeeze(Momentsest(20,1,index1,index2,[vss_set])))./squeeze(Momentsest(18,1,index1,index2,[vss_set]))];
buildupasym(4,:,3)=[squeeze(Momentsest(10,1,index1,index2,[vss_set]))];

subplot(2,1,find(delta_calib==delta_calibvec))
 hold on
 plot(xx'+1, 1-beta_dcdy_data*ones(length(xx),1), ':','LineWidth',2.5,'Color','k');
plot(xx'+1,1-buildupasym(4,:,2)','linestyle', '--','LineWidth',2.5,'Color','k');
if plotasym==1
yyaxis right
plot(xx'+1, buildupasym(3,:,2)','linestyle', '-','LineWidth',2.5,'Color','k');
end
end

subplot(2,1,1)
title('Moderate insurance','Interpreter','Latex','FontWeight','bold','FontSize',24)%['$\delta=$' num2str(delta_calibvec(1))],'Interpreter','Latex','FontWeight','bold','FontSize',28)
 yyaxis left
ylim([0.4,1])
 yyaxis right
ylim([0,0.37])
xlim([2,vss_max])

subplot(2,1,2)
title('Strong insurance','Interpreter','Latex','FontWeight','bold','FontSize',24) %['$\delta=$' num2str(delta_calibvec(2))],'Interpreter','Latex','FontWeight','bold','FontSize',28)
if plotasym~=1
h1= legend('Degree of insurance, data','Degree of insurance, model','Location','Best','FontWeight','bold')%,'$Var_{dc|dy>0}/Var_{dc|dy<0}$ (rhs)'  
else
h1= legend('Degree of insurance, data (lhs)','Degree of insurance, model (lhs)','Asymmetry of cons response (rhs)','Location','Best','FontWeight','bold')%,'$Var_{dc|dy>0}/Var_{dc|dy<0}$ (rhs)'
end
set(h1,'Interpreter','Latex','FontSize',16.5,'Location','SouthOutside','Orientation','Horizontal')
 yyaxis left
ylim([0.4,1])
 yyaxis right
ylim([0,0.37])
xlim([2,vss_max])
xlabel('Village Size','FontSize',24,'FontWeight','bold','Interpreter','Latex')

set(figure2,'Units','Inches');
pos = get(figure2,'Position');

set(figure2,'PaperPositionMode','Auto','PaperUnits','Inches','PaperSize',[pos(3), pos(4)])
set(gcf, 'PaperOrientation', 'landscape')
if write==1
print(figure2,[currentfolder '\Latex\Moments_fct_of_Vil_size_cond_' date '_paper'],'-dpdf','-r0')
end
pause

%  ================================
%  Figure 3
%  ================================
clf
close all

figure3=figure
subplot(1,3,1)
title('CD model','Interpreter','Latex','FontWeight','bold','FontSize',24)%['$\delta=$' num2str(delta_calibvec(1))],'Interpreter','Latex','FontWeight','bold','FontSize',28)
load([currentfolder1 '\Scattersome_name1110010.mat'])
yy_cond_plot=yy_cond(:,1:2000);
cc_cond_plot=cc_cond(:,1:2000);
plot(yy_cond_plot(:),cc_cond_plot(:),'ko')
xlabel('Income, log-difference','FontSize',16,'FontWeight','bold','Interpreter','Latex')
ylabel('Consumption, log-difference','FontSize',18,'FontWeight','bold','Interpreter','Latex')

subplot(1,3,2)
title('ID model','Interpreter','Latex','FontWeight','bold','FontSize',24)%['$\delta=$' num2str(delta_calibvec(1))],'Interpreter','Latex','FontWeight','bold','FontSize',28)
load([currentfolder1 '\Scattersome_name1010010.mat'])
yy_cond_plot=yy_cond(:,1:2000);
cc_cond_plot=cc_cond(:,1:2000);
plot(yy_cond_plot(:),cc_cond_plot(:),'ko')
xlabel('Income, log-difference','FontSize',16,'FontWeight','bold','Interpreter','Latex')
ylabel('')


subplot(1,3,3)
title('SI model','Interpreter','Latex','FontWeight','bold','FontSize',24)%['$\delta=$' num2str(delta_calibvec(1))],'Interpreter','Latex','FontWeight','bold','FontSize',28)
load([currentfolder1 '\Scattersome_name1210010.mat'])
yy_cond_plot=yy_cond(:,1:2000);
cc_cond_plot=cc_cond(:,1:2000);
plot(yy_cond_plot(:),cc_cond_plot(:),'ko')
xlabel('Income, log-difference','FontSize',16,'FontWeight','bold','Interpreter','Latex')
set(gcf,'PaperPositionMode','Auto','PaperUnits','Inches','PaperSize',[pos(3), pos(4)])
set(gcf, 'PaperOrientation', 'portrait')
if write==1
print(gcf,'-loose',[currentfolder '\Latex\scatterGE_IDvsCDvsSI_' date],'-depsc2','-r0')
end


% Figure 1
% replicate graph with rho and beta as in LTW
 figure(1)
title('Consumption and income growth under limited commitment','Interpreter','Latex','FontWeight','bold','FontSize',26)%['$\delta=$' num2str(delta_calibvec(1))],'Interpreter','Latex','FontWeight','bold','FontSize',28)
 load([currentfolder '\ScatterVilcoal_15_5_19_vill_cd_incproc1010010.mat'])
 yy_cond_plot=yy_cond(:,1:2000);
 cc_cond_plot=cc_cond(:,1:2000);
 plot(yy_cond_plot(:),cc_cond_plot(:),'ko')
 xlabel('Income, log-difference','FontSize',24,'FontWeight','bold','Interpreter','Latex')
 ylabel('Consumption, log-difference','FontSize',24,'FontWeight','bold','Interpreter','Latex')
